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Resumo 


No presente trabalho faz-se uma abordagem da aplicação do Cálculo Fracio¬ 
nário na modelagem de sistemas físicos complexos. Começando pelo processo de 
difusão de soluções e terminando com a modelagem do processo de difusão de nêu¬ 
trons em reatores nucleares. Para isto, necessita-se fazer uma abordagem histórica 
e conceituai do Cálculo de Ordem Fracionária, assunto que requer certo grau de 
domínio conceituai especíüco e refino matemático. Para estudo inicial deste, se faz 
necessário o entendimento de algumas funções especiais, sendo elas: Gama, Beta 
, Mittag-Leffler, Função R e Função F. Dentro do Cálculo Fracionário existem di¬ 
versas formulações, neste texto serão estudadas as de Caputo, e Riemann-Liouville, 
sendo essas as mais conhecidas e utilizadas. Dentro dessas formulações será feita a 
diferenciação entre integral e derivada fracionária para facilitar o estudo. 

Palavras chave: Cálculo Fracionário, difusão anômala, nêutrons, difusão de nêu¬ 
trons, Reatores nucleares, Gama, Beta, Mittag-Leffler, Caputo, , Riemann-Liouville. 

Abstract 

In the present work an approach is applied to the application of the Fractional 
Calculus in the modeling of complex physical Systems. Starting with the process 
of diffusion of Solutions and ending with the modeling of the process of diffusion 
of nêutrons in nuclear reactors. For this, it is necessary to make a historical and 
conceptual approach of the Fractional Order Calculus, subject that requires certain 
degree of specihc conceptual domain and mathematical rehnement. For the initial 
study of this, it is necessary to understand some special functions, being: Ganima, 
Beta, Mittag-Leffler, Function R and Function F. Within the Fractional Calculus 
there are several formulations, in this text will be studied those of Caputo, and 
Riemann -Liouville, these being the most known and used. Within these formula¬ 
tions will be made the differentiation between integral and fractional derivative to 
facilitate the study. 

Key words: Fractional Calculus, Anomalous Diffusion, Nêutrons, Nêutron Diffu¬ 
sion, Nuclear Reactors, Gamma, Beta, Mittag-Leffler, Caputo,, Riemann-Liouville. 
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Capítulo 1 
Introdução 


1.1 Difusão 

A difusão pode ser definida como o movimento randômico de partículas de re¬ 
giões de maiores concentrações para regiões de menores concentrações de maneira a 
diminuir o gradiente de concentração de partículas no meio. 


1.2 Contexto Histórico 

As primeiras medidas sistemáticas de um processo de difusão foram feitas pelo 
químico escocês Thomas Graham (1805-1869), que estudou difusão de gases e líqui¬ 
dos entre 1828 e 1850. Por exemplo, para estudar a difusão de substâncias na água, 
Graham encheu dois reservatórios com diferentes concentrações de ácidos ou sais em 
água, permitiu que os conteúdos dos dois reservatórios se comunicassem por vários 
dias, e então analisou as respectivas concentrações. 

Graham observou que a quantidade de soluto transportada decrescia em fun¬ 
ção do tempo e que a quantidade de soluto que se difundia durante um período 
üxo de tempo de um reci-piente cheio para um recipiente vazio era proporcional à 
quantidade de soluto inicialmente presente. 

Em 1855, o médico e üsiologista alemão Adolf Fick (1829-1901), então pro¬ 
fessor associado de anatomia e fisiologia da Universidade de Zurique, propôs leis 
fenomenológicas para a difusão de um soluto em um solvente (que atualmente le¬ 
vam o seu nome) análogas aos modelos de Fourier para a condução de calor e de 
Ohm para a condução de eletricidade (mostrando que um tipo de processo pode ser 
comum a vários sistemas físicos).[1] 
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1.3 DIFUSÃO EM ESTADO ESTACIONÁRIO 


O fluxo de partículas (p(x, t ) é número líquido de moles de partículas cruzando 
uma área unitária perpendicular ao eixo x, localizada em x, por unidade de tempo, 
no instante t. As unidades de 4>(x,t) são mol/(cm 2 .s). 

Em três dimensões, o fluxo é definido como o número líquido de moles de 
partículas que cruza uma superfície de área unitária em uma unidade de tempo. 
Portanto, o fluxo é um vetor (isto é, ele tem componentes nas direções x, y e z). Em 
uma dimensão, o fluxo reduz-se a um escalar. 

Quando o fluxo (j)[x , t) não varia com o tempo (a concentração C também 
não varia com o tempo) e temos a DIFUSÃO EM ESTADO ESTACIONÁRIO (ou 
DIFUSÃO EM REGIME PERMANENTE). 

Para que fluxo (cp(x,t)) não varie com o tempo é necessário que (ft(x,t) 
também não varie com a posição. Assim, 

dC 

cp(x, t ) = cte. —» D—~ = cte. —> C(x) = ax + ò, a > 0 (1.1) 

ox 

1.3.1 A Primeira Lei de Fick 

Com base no estudo sobre o fluxo de calor entre dois corpos a temperaturas 
diferentes, Fick propôs que o fluxo por difusão (p(x , t), do soluto S é proporcional ao 
gradiente de concentração entre os dois lados da membrana semipermeável. Em que 
defini-se como membrana semipermeável aquelas membranas que deixam difundir 
apenas o solvente, impedindo a difusão do soluto. A imagem 1.1, abaixo, ilustra a 
situação formulada. 
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Membrana semipermeável (S) 


Figura 1.1: Fluxo de particulas através de uma membrana semipermeável[1] 

O fluxo por difusão é expresso na equação (1.1), na qual C(x,t) é concentração 
em função da posição(x) e do tempo (t), D é a constante de difusão, x é variável 
que representa o fluxo em uma dimensão e t representa o tempo : 
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<f)( x ,t) =-D—, (1.2) 

O sinal negativo na equação acima indica que o fluxo ocorre no sentido contrário 
ao do gradiente de concentração, isto é, no sentido do meio de maior concentrações 
para o de menor concentrações. 

Contudo, é importante lembrar, que a Lei de Fick é apenas uma das muitas 
formas da equação do transporte. A figura 1.2, abaixo, apresenta as principais 
variações da equação de transporte: 


Various fornis of the transport equation. 
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Figura 1.2: Equações de transporte [2] 

A equação diferencial do processo de difusão implica que as partículas obe¬ 
decem a derivada de ordem um em relação ao tempo. Os processos que fogem dessa 
relação são conhecidos como anômalos, esse conjunto de processos difusivos anôma¬ 
los é dividido em dois: processos subdifusivos e super difusivos, que serão abordados 
a frente com maior profundidade. [2] 
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1.4 DIFUSÃO EM ESTADO NAO-ESTACIONARIO 


Em sua maior parte, as situações práticas, envolvendo difusão ocorrem em condi¬ 
ções de ESTADO NÃO-ESTACIONÁRIO (ou REGIME TRANSITÓRIO ou CON¬ 
DIÇÕES TRANSIENTES). 

No processo de difusão em estado não-estacionário, tanto o fluxo de difusão, 
quanto o gradiente de concentração, numa dada posição x, variam com o tempo t. 
Como consequência, ocorre um ajuntamento ou esgotamento líquido do componente 
que se encontra em difusão. 

1.4.1 A relação de Einstein entre Difusão e Viscosidade 

Einstein adotou desde cedo uma visão realista, objetiva, sobre a existência de 
átomos e moléculas. Na sua tese de doutoramento Einstein analisa o fenômeno 
de difusão das partículas do soluto numa solução diluída (moléculas de açúcar em 
água) com o objetivo de obter estimativas para o número de Avogadro e o diâmetro 
das partículas do soluto. As propriedades termodinâmicas das soluções diluídas já 
tinham sido suficientemente estabelecidas. Sabia-se, por exemplo, que a pressão 
osmótica, exercida pela solução sobre uma membrana semi-permeável, impedindo 
a passagem do soluto, comporta-se de acordo com a lei dos gases perfeitos. Na 
parte inicial da tese, Einstein faz um cálculo hidrodinâmico, com base nas equações 
de Navier-Stokes para o escoamento de um fluido incompressível, a hm de obter a 
viscosidade efetiva do fluido na presença do soluto. No modelo adotado, as moléculas 
do soluto são esferas rígidas, não interagentes, e bem maiores do que as moléculas do 
solvente. O resultado final, que ainda tem um pequeno erro e precisou ser corrigido 
alguns anos depois por um aluno de Einstein, é dado por[4] 

rf* = 7/(1 + $) (1.3) 

onde rj* é a viscosidade efetiva, rj é a viscosidade do solvente puro, e $ é a fração 
do volume total ocupado pelas partículas do soluto. Supondo que o soluto seja 
constituído por partículas esféricas de raio a, a fração de volume é dada por, 

$ = (1.4) 

onde Ns é o número de moléculas do soluto e V é o volume total da solução. Uti¬ 
lizando a densidade de massa, p, e a massa molar do soluto, m, que são grandezas 
experimentalmente acessíveis, temos 

pV 

N s = '—N a (1.5) 

m 
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onde Na é o número de Avogadro. A partir dessas expressões, podemos escrever a 
relação 

^-7ra 3 N s = — - 1 ( 1 . 6 ) 

óV p 

Já qne as viscosidades podem ser medidas, qne a densidade de massa e a massa 
molar do soluto são conhecidas, aparecem como incógnitas o raio a das partículas 
do soluto e o número de Avogadro Na- Na parte final da tese, Einstein recorre a um 
argumento engenhoso, deduzido de forma alternativa no artigo sobre o movimento 
browniano, a fim de obter uma segunda relação entre a e Na- Vamos discutir esse 
argumento de Einstein. A idéia engenhosa consiste em considerar uma força K (na 
direção do eixo x) atuando sobre as partículas grandes da solução (as moléculas de 
açúcar são muito maiores do que as moléculas de água), contidas em um volume 
elementar de comprimento Are e seção transversal As. Para esferas rígidas de raio a 
e velocidade v, mergulhadas num fluido com viscosidade rj, essa força deve ser dada 
pela lei (de atrito viscoso) de Stokes, 

K = Girarjv (1.7) 


de acordo com o texto de mecânica de Kirchhoff, que Einstein estudou durante os 
seus anos de formação em Zurique (e que é a única referência citada no primeiro 
trabalho sobre a teoria do movimento browniano!) Mas as partículas se difundem 
pelo fluido devido ao gradiente de pressão. Então, nesse volume elementar, estariam 
sujeitas a uma força por unidade de volume, ao longo do mesmo eixo x, dada por 

d P P(x + Ax) — P(x ) 

= fim -1- 

OX Az-s>0 /\x 


Podemos agora escrever uma equaçao de balanço entre essas duas forças, 

m d P 


K = 


= Gnarjv 


(1.9) 


pNa dx 

onde pé a densidade de massa e m é a massa molar do soluto. Obtemos assim 
uma expressão para a velocidade v das partículas, que nos remete ao fluxo (f) ao 
longo do eixo x (quantidade de massa das partículas atravessando a seção de área 
AS durante o intervalo de tempo Aí), 


d> = pv — 


m d P 
pN A dx 


= Girarjv 


( 1 . 10 ) 


Pela equação geral dos gases ideais para pressão osmótica, sabendo que ^ — p, 


temos 


PV = nRT 


P = 


RTp 


m 


( 1 . 11 ) 



que Einstein justificava através de argumentos de física estatística para um sistema 
de partículas não interagentes, temos 


<f> = pv 


m d P 
pN A dx 


Gna/qv 



( 1 . 12 ) 


de onde vem 


D 


RT 

QnpaN A 


(1.13) 


que é uma das expressões conhecidas de Einstein, relacionando o coeficiente de difu¬ 
são D, experimentalmente acessível, com a temperatura e a viscosidade do fluido.[4] 


1.4.2 A segunda Lei de Fick 

Para modelar o fenômeno da difusão em estado não-estacionário, é utilizada a 
equação diferencial parcial (EDP) 

dC{x,t) _ d í dC(x,t) \ 

dt dx V dx ) 1 j 


conhecida por Segunda lei de Fick 

Se o coeficiente de difusão não depende da composição (portanto, da posição), 
a segunda lei de Fick se simplifica para 

dC(x,t) d 2 C(x,t ) 

—dt r- = (L15) 

Principais condições para aplicação da segunda Lei de Fick 
Para o caso hipotético, em que esteja sendo estudado, um certo volume de uma 
substância de concetração tal, que atenda as condições abaixo: 

1. O perfil de concentração muda com o tempo: o gradiente de concentração não 
é constante no interior do volume considerado. 

2. Quando: <j>(x,t) Entrada > Saída há acúmulo da espécie difundida no vo¬ 

lume considerado. 
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1.5 O Cálculo Fracionário 


1.5.1 Histórico 

Atualmente, tem se observado um grande interesse, por parte da comunidade 
acadêmica dos cientistas e engenheiros, pela modelagem matemática de fenômenos 
com o chamado Cálculo de Ordem não Inteira ou popularmente conhecido como 
Cálculo Fracionai ou Cálculo Fracionário (CF). 

Para muitos estudiosos o (CF) nasceu na troca de carta datada de 30 de 
setembro de 1695, entre o Marquês de L’Hôpital( 1661-1704) e Leibniz (1646-1716) 
na qual o assunto derivada foi considerado de forma bastante interessante, a saber: 
Leibniz formulou uma questão envolvendo uma possível generalização da derivada 
de ordem inteira para uma ordem, em princípio, arbitrária. L’Hôpital devolveu 
a pergunta para Leibniz, questionando-o no caso particular em que a ordem da 
derivada (equação 1.16) fosse meio (n = |). 


d n f(x) 

dx n 


(1.16) 


Em uma audaciosa e profética resposta, Leibniz apresenta o resultado e 
afirma: isto é, aparentemente, um paradoxo que um dia vai gerar várias consequên¬ 
cias importantes. [5] 


1.5.2 Funções Especiais 
A Função Gama 

A Função Gama, ou função de Euler do Segundo Espécie, é uma extensão do 
significado do fatorial de um número natural para o fatorial de um número fracio¬ 
nário, real, ou mesmo complexo. 

Para melhor conhecer suas propriedades, vamos defini-la por três formas diferen¬ 
tes: 


Primeira definição - Através do limite infinito (Euler): 

1 • 2 • 3...n 


r(z) = lim 


-n 


z #0,-1, -2, -3... 


(1.17) 


n-s-oo z(z + l)(z + 2 )...(z + n) 

Segunda definição - Através da integral definida (também devida a 
Euler): 


r(z) = 


e H z 1 dt , 


3fi(z) > 0 


(1.18) 
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Propriedades 


r (n + 1 ) = n\. 

(1.19) 

T{z + l) = zT{z) 

( 1 . 20 ) 


Função Beta 

A função Beta ou função de Euler de Primeira Espécie apresenta significativa 
relevância no Cálculo Fracionário. A definição na forma integral é: 

/3(p,q) = í t pL (l — t) q ~ 1 dt, > 0 , 3Ç(g) > 0 . (1-21) 

J o 

Agora a forma trigonométrica, obtida quando se substitui a variável de integração 
t por cos 2 9 na equação 2.17, que será muito útil em demonstrações posteriores. 

7r 

/3(p,q) = 2 í cos 2p ~ 1 9sen 2q ~ 1 9d9, 5í(p) > 0 , !R(g) > 0 . (1-22) 

J o 

Também se conhece a forma racional quando se substitui a variável de integração 
t por na equação 2.17: 


P(p, q) 



u p - 1 

7 -- ——du, 

(1 + u)p+i 


3?(p) > 0 , 9ff(ç) > 0 . 


(1.23) 


Outra forma importante é a generalização da função Beta, em que o limite su¬ 
perior passa a ser um número arbitrário a, pode-se obtê-la a partir da subtituição 
da variável de integração t por 

y p - l (a - y) q ~ l dy, K(p) > 0 , &(g) > 0. (1.24) 

determinada relação com a função Gama: 


A função Beta possui 


P(p,q) 


r(p)r(g) 
r ip + q) ' 


(1.25) 
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1.5.3 A Integral Fracionária 

Como demonstrado na introdução, a conceituação sobre integrais e derivadas de 
ordem fracionária foi uma consequência natural dessas com ordem inteira. Indaga¬ 
ções a respeito, por exemplo, da derivada de ordem meia de uma função, estimularam 
o estudo do cálculo fracionário. 

A comum formulação para a integral fracionária pode ser obtida diretamente 
através de repetidas integrações de uma função. Esta aproximação é conhecida como 
abordagem de Riemann- Liouville. Abaixo se verifica a fórmula para n-integrações 
de uma função f(t). 

í ■[ f( T )dr = 7 — vry f (t - t)^ 1 f {r)dr (1.26) 

J t n —\ J o o 

A fim de abreviar a fórmula, é introduzido o operador J n como em: 

= -An A* - (1.27) 

O operador J n pode ser trocado por mas é preferível usar o primeiro de 

modo a evitar possíveis conflitos, uma vez que se usa o operador D n nas derivadas 
fracionárias. 

A equação (1.23)’e limitada, pois n se restringe a ser um inteiro. A principal 
limitação se dá na utilização do fatorial, que não tem significado para valores não- 
inteiros. Tentando estender a aplicabilidade desta equação, lança-se mão da função 
gama no lugar do fatorial. Assim, é possível generalizar para todo aeR +1 


J n f(t) 


1 

FR 


(t — t ) 11 - 1 f {t)(Lt 


(1.28) 


Propriedades 

Algumas propriedades são observadas nas integrais fracionárias. Em primeiro 
lugar, considera-se que a integração de ordem a = 0 é o operador identidade. 

Df(t) = f(t ) (1.29) 

Uma segunda propriedade, baseada na definição natural de integral e no princípio 
do qual ela veio ( Cauchy ), é possível dizer que: 

n J m = J n+m = J m J n , m e neR. (1.30) 

então, 

jajp = ja+p = jPja ? ^ £ R (L 31 ) 
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A linearidade também é uma das características dessas integrais fracionárias. 

J a [afi(t ) + bf 2 (t)\ = aJ Q h{t ) + bJ a f 2 (t ) (1.32) 

Uma propriedade adicional aparece após a introdução da função f I> 0 . 

=> ®a(t) * f(t ) = í — ^ — f(r)dr (1.33) 

T(a) J o r(a) 

Assim, tem-se: 

J a f(t) = $<*(*) * /(í) = ^ 7 — Í (t ~ T) a - l f(j)dT (1.34) 

r(a) J o 

Nota-se através da equação anterior que a integral fracionária pode ser expressa 
como a convolução de dois termos, <& a e f(t). A partir da seguinte transformada de 
Laplace: 


£{<"} = *,«>0 ( 1 . 35 ) 

é possível obter que: 

£{*“-’} = ^ ( L36 ) 

e também: 

£{$(«)} = 4 (1.37) 

1.5.4 A Derivada Fracionária 

A partir da formulação de Riemann-Liouville para a integração de ordem arbi¬ 
trária, pode- se estabelecer uma definição para a derivação de ordem fracionária. A 
denominada derivada de ordem fracionária de Riemann-Liouville está relacionada à 
operação inversa da integração: 


Assim, 


D a m = D m [rf{t)} , t> 0 


(1.38) 


D a m 


d m f 1 [t f{r) J ] 
dt m Lr(m-a) J o (t—r )“+ 1 m J 


d m 

dt m 


m 


m — 1 < a < m 

(1.39) 

a = m 
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A explicação do procedimento é simples. Tendo encontrado o inteiro m, em 
que m — 1 < a < m, o primeiro passo é integrar a função f(t) por ordem rri — a 
e, em seguida, diferenciar a função resultante f(t) pela ordem m. Desta maneira, 
obtendo-se uma diferenciação de ordem resultante a. A equação (17) estabelece a 
definição de D a f{t) e é conhecida como definição do lado esquerdo: LHD-Left Hand 
Definition. 

Existe uma segunda definição que é conhecida como definição do lado direito. 
Esta usa as mesmas operações que a anterior, porém em ordem inversa. A RHD- 
Right Hand Definition é dada por: 


D a J{t) = 


1 [* f m {r) 

r (m - a) J 0 (t - T ) a+1 - ,r 
d m 

— — f(t), a — m 
dt m v ’ 


-dr, m — 1 < a < m 


(1.40) 

(1.41) 


Esta última, apesar de ser referida aqui como RHD, foi originalmente formu¬ 
lada por Caputo, e é por isso, com um ente conhecida como derivação fracionária de 
Caputo. 


1.5.5 Caputo x Riemann-Liouville 

Inicialmente, foi suposto uma função f(t) na qual foi aplicada a derivada fraci¬ 
onária de ordem a, seguindo as duas definições anteriores: D a {t) e Df(t). Estas 
definições segundo Riemann-Liouville e Caputo, respectivamente, só serão equiva¬ 
lentes se f(t) for uma função causal, isto é, é nula para t < 0. 

Outro aspecto que deve ser observado é que embora para polinómios não-constantes 
as definições segundo Riemman-Liouville e Caputo coincidam, isto não ocorre para 
constantes, uma vez que na definição de Caputo primeiramente toma-se a derivada 
de ordem inteira e depois a integral fracionária. A derivada de ordem fracionária de 
uma constante segundo Caputo é nula, ou seja: Dfl = 0, o que não ocorre quando 
se calcula a derivada a partir de Riemann-Liouville. Esta última fornece: 

D ° c = yfif—— \ ( L42 ) 

1 (1 — a) 

Considerando a metodologia da transformada de Laplace nas duas definições, é pos¬ 
sível fazer algumas observações. 

Riemann-Liouville: 
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Caputo: 


m =1 

C{D a f{t)} = s°F(s ) - ^ s k D ía - k ~ 1) f{ 0) (1.43) 

k=0 


m =1 

£{D“/(f)} = s Q F(s) - ^ s (a-fc- 1 )^(fc)( 0 ) > «eC, m E N (1.44) 

k=0 

A transformada de Laplace segundo Caputo parece ser mais apropriada na reso¬ 
lução de equações diferencias fracionárias, pois requer o conhecimento de condições 
iniciais da função e suas derivadas de ordem inteira, que têm interpretações físicas. 

Neste contexto, embora a definição de Riemann-Liouville para derivadas fra¬ 
cionárias seja mais conhecida, a formulação de Caputo tem maior aplicabilidade, 
mesmo sendo esta última mais restritiva. 


1.5.6 Função de Mittag-Leffler 

Juntamente com as funções anteriores definidas por Euler, a Função de Mittag- 
Leffler (na verdade, toda uma classe de funções), é de presença ubíqua no CF. A 
Função de Mittag-Leffler trata-se de uma generalização da exponencial, e tem um 
papel equivalente na resolução de Equações Diferenciais Fracionárias lineares ao que 
a função exponencial tem para as Equações Diferenciais lineares de ordens inteiras. 
A Função de Mittag-Leffler foi introduzida primeiramente em 1903 pelo matemático 
sueco Magnus Gosta Mittag-Leffler, e posteriormente generalizada contendo mais 
parâmetros. Essa função é dependente de um parâmetro complexo a, como se 
mostra na série a seguir [6]: 


E a (x) — 1 + ———— + ———— + ———— + ... + zvr-—T + ••• (1 

r(a -|- 1) r(2a -f- 1) r(3a -|- 1) L(ttq; -|- 1) 

Para o caso particular onde a = 1, a função de Mittag-Leffler se reduz à própria 
exponencial. 



Ei(x) = J2 


x 


-í r(n + 1) 

n=0 v 7 


= e 


(1.46) 


Essa função foi generalizada para dois parâmetros, originalmente por Agarwal 
da seguinte maneira: 


E a A x ) = 


X 


n^r (a.n + 0) 


Sft(of) > 0,9ft(/3) > 0. 


(1.47) 
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1.5.7 A Função R e a Função F 

Desde 1903, a função Mittag-Leffler foi modificada para ter uma representação geral 
para a resolução de equações diferenciais fracionárias. Ultimamente, a tentativa de 
Miller e Ross (1993), Lorenzo e Hartley (1999) para generalizar ainda mais resultou 
em duas funções importantes: a função R e a função F, listadas abaixo[7]: 


a n t( n+ V q - 1 

,M “ 5 r({ n + i} 9 ) 

(1.48) 

£{*>,«]} = s „_ a 

(1.49) 

r 1 ^ a n t (n+l)q-l-v 

~ “o r W”+m— v) 

(1.50) 

£{R q , v [a,t]} = 

si — a 

(1.51) 

1.6 Equação da difusão anômala modelada via equa¬ 
ção da difusão fracionária 

Difusão unidimensional em um meio finito [6] 

A difusão convencional é representada pela equação (1.51) . 
fracionária: 

d a C(x,t ) d 2 C(x,t) 

dt a ~ dx 2 ’ 

, no caso da difusão 

(1.52) 


onde C(x,t) representa a concentração de solnto no ponto x e no tempo t. 

Pelos estudos das equações diferenciais parciais unidimensionais, se sabe que, 
em geral, para a resolução destas são necessárias duas condições de contorno e uma 
inicial. O problema tratado aqui será referente a difusão em um meio unidimen¬ 
sional, assim a função desejada irá variar apenas com relação ao eixo x e com o 
tempo. O objeto de estudo terá comprimento total igual a L, onde sua matéria está 
totalmente distribuída de acordo com uma função genérica c(x) no instante inicial[6]. 
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Essas condiçoes de contorno mostram que o meio é fechado, portanto não há 
fluxo nas extremidades. 

Se faz a hipótese típica do método de separação de variáveis, na qual a função 
desejada pode ser representada como um produto de duas funções, na qual uma 
depende apenas de x e a outra depende apenas de t: 

C(x,t) = X(x)T(t). (1.54) 

Substituindo a hipótese na equação diferencial de difusão fracionária, tem-se: 

X(x)T {a \t) = DX"(x)T(t). (1.55) 

Pelo resultado acima, pode-se perceber que o método de separação de variáveis 
funcionou, onde se consegue isolar o que depende da variável x de um lado e o 
que depende de t do outro, assim tudo isso se torna igual a uma constante, pois 
é a única solução possível para uma igualdade entre dois membros com variáveis 
independentes de modo que possa ser verdade para qualquer alteração entre estas, 
como se pode ver abaixo: 

Thd(í) X"(x) , , . . 

mi) (consi “ nfe )' (L56) 

Agora sim, o problema que era a resolução de uma EDP se torna a resolução de 
duas EDO’s. Nesse momento tem-se três hipóteses para a constante —À, esta pode 
ser negativa, positiva ou nula. Estudando o caso em que À é igual a zero, tem-se: 

X"(x) = 0. (1.57) 

A solução mais geral para a X(x) nesse caso será: 

X(x) = Bx + C, sendo B e C constantes quaisquer. (1.58) 

Utilizando a primeira condição de contorno, equação (1.53): 

X'(0)T(t) = 0. (1.59) 

Para isso ser verdade e ser considerado solução, X’(0) tem que ser igual a 0, 
pois, caso T(t) seja 0, a equação da concentração será sempre nula, como pode 
ser facilmente veriücado substituindo essa condição na equação (1.54), o que só 
seria solução do problema para o caso particular e trivial onde c(x) = 0, o que não 
interessa, pois, se procura uma solução geral. Aplicando X’(0)=0 em X(x), tem-se 
que [6]: 
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logo: 


X'(0) = B = 0, 


(1.60) 


X(x) = C. (1.61) 

Utilizando a segunda condição de contorno: 

X\L)T{t) = 0. (1.62) 

Como explicado anteriormente só seria solução se X’(L) for 0, o que leva a: 

X(L) =C. (1.63) 

Assim se chega ao resultado de que a concentração não varia com a posição, o 
que não é interessante. 

Estudando os outros casos, —À > 0 e — À < 0, tem-se as seguintes EDO’s: 

T (q) (í) = -A DT(t), (1.64) 

X"(x) = -\X(x). (1.65) 

Para os dois casos tem-se que a equação relacionada a x pode ser resolvido a 
partir de um palpite inicial (ansatz), este palpite será uma exponencial em razão 
do aparecimento constante desta função na resolução de equações diferenciais or¬ 
dinárias, assim supõe-se que X(x) será do tipo de rx , sendo der constantes. Para 
encontrar os valores de r, substitui-se o palpite na equação (1.65), obtendo: 

r 2 de rx = -A de rx , com d G & - (0). (1.66) 

Dividindo toda equação por de rx , que pode-se ter certeza que não será um termo 
nulo, tem-se: 

r 2 = -A, (1.67) 

logo: 

r = ix/^Ã. (1.68) 

Considerando que — A > 0, tem-se que: 

X(x) = c 1 e V ^ Xx + c 2 e-^ Xx . (1.69) 

Substituindo a primeira condição de contorno, X’(0) = 0 utilizada anterior mente, 
percebe-se que só é verdade caso d\ = <á 2 - Assim: 

X{x) = + e _>/ZSa: ]. (1.70) 
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Com a segunda condição de contorno, X’(L) = 0, obtêm-se que há solução ape¬ 
nas se cl\ for 0, que não interessa, pois leva novamente ao caso trivial. 

Agora resta o último caso, no qual —À < 0. Para este pode-se partir da equação 

(1.74) , na qual, fazendo o uso de números complexos: 

X(x) = die iV * x + d 2 e~ iVXx . ( 1 . 71 ) 

Assim, como se pode perceber, essa equação é uma soma de um seno e um 
cosseno. 

X(x) = diCos(\f\x) + d 2 sen(VXx). (1-72) 

Utilizando a primeira condição de contorno X’(0) = 0, encontra-se que d 2 é igual 
a zero. Logo: 

X(x) = d 1 cos(VXx ). (1-73) 

Utilizando a segunda condição de contorno X’(L) = 0, encontra-se que a f\L deve 
ser múltiplo inteiro de n: 

'Yl^TT^ 

X = , com n G N. (1-74) 

Como n pode ser qualquer inteiro, esse problema possui infinitas soluções, sendo 
que a cada uma delas se associa uma nova constante. Substituindo o À da equação 

(1.74) e considerando que n pode ser um número natural qualquer, tem-se que a 
função X(x) é igual a um somatório de cossenos na sua forma mais geral, pois, uma 
combinação linear de soluções também é solução. 

OO 

X(x) = y^ j d n cos( 1 ^-'] (1.75) 

n= 0 

Agora que foi encontrada a função relacionada a x, serão feitos os cálculos para 
t a partir da equação (1.64). Para a resoluçõo da EDO relacionada ao tempo, 
aplica-se a transformada de Laplace nos dois membros, em razão de À e D serem 
constantes, estes não sofrem transformação. 

C{T°(t)} = -X DC{T(t)}. (1.76) 

Para um estudo bem genérico, será necessário neste momento a distinção entre 
dois processos anômalos, os subdifusivos e os super difusivos. Os processos subdifu- 
sivos são aqueles em que a ordem da derivada varia de 0 até 1, já os superdifusivos 
são aqueles em que a ordem da derivada varia de 1 até 2, essa diferença promove 
mudanças no momento da resolução. 
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Quando se aplica a transformada de Laplace na parte esquerda da equação (1.64) 
de um processo subdifusivo, faz se necessário o uso de uma forma genérica da regra 
da transformada de Laplace da derivada pela formulação de Caputo. Para o caso 
particular, tem-se que: 


C{T a (t)} = s a C{T(t)} -s“ _1 T(0). (1.77) 

Substituindo na equação (1.76) e isolando a transformada de Laplace de T(t), 
obtêm-se: 

C{T(t)} = --(1.78) 

L WJ s“ + DA v ’ 

Agora aplica-se a transformada de Laplace inversa na equação anterior, mas para 
encontrar um resultado que pode passar por uma melhor análise será utilizado um 
segundo resultado expresso no apêndice deste trabalho, referente a transformada de 
Laplace de um caso geral da função de Mittag-Leffler. Fazendo uso deste, obtêm-se 
que: 

T(t) = T(0)E a i(—DXt a ). (1.79) 

Substituindo o valor de À e considerando n qualquer número natural, pode-se 
escrever a solução anterior como uma combinação linear de soluções: 

m = • (1.80) 


Assim se obtêm a solução geral do problema fazendo a união da solução X(x) e 

T(t): 


C{x,t) = 


= 2 ^ 9nCOS 
n =0 


/ rmx \ 

\~L 


j E a \ 


-Dn 2 n 2 t a 
L 2 


(1.81) 


Quando se aplica a transformada de Laplace na parte esquerda da equação (1.76) 
de um processo superdifusivo, o resultado será o seguinte: 


C{T a (t)} = s°C{T(t)} — s a ~ 1 T(0) — s a ~ 2 T'(0). 


(1.82) 


Substituindo na equação (1.76) e isolando a transformada de Laplace de T(t), 
obtêm-se: 

c[rm = >_ A3 + 

1 J s" + DX s a + D X 

Com o auxílio da transformada inversa e dos resultados do apêndice, se consegue 
a seguinte solução: 


T(t) = f n T(0)E aA 


-Dn 2 7T 2 t a 

L 2 


+ fnT\ti)tE a 2 


-Dn 2 n 2 t a 
L 2 


(1.84) 
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Agora que já se tem os resultados desejados, se consegue a solução geral para a 
concentração: 



(1.85) 


Utilizando agora a condição inicial, C(x,0) = c(x), pode-se encontrar os g n ’s a 
partir dos conhecimentos das séries de Fonrier 1 . 

c(x) = C(x, 0) = £ 9nCos(——jr ^-). (1.86) 

n=0 

Desse modo: 

1 f L n%x 2 f L nirx 

9n=y c(x)cos(—-—)dx = — / c(x)cos( —-— )dx, n — 0,1,2,3... (1-87) 

L J-L L L J o L 

Essa sim é uma resposta bem geral para este problema, sendo qne a função da 
condição inicial irá sofrer uma extensão par em razão de ser composta apenas por 
cossenos. 


A constante K qne é aqui empregada é uma constante real igual a a variável 
t inicia-se em 0, a constante C se refere a concentração inicial (T(0)) e a constante D 
se refere a variação da concentração inicialmente (T’(0)). Para a subdifnsão, nossa 
concentração será analisada em x=0, aqui representada por y: 


y(t) = C Y1 

n=0 


_ \)n K 2n 

r {na + 1) 


0 < a < 1. 


Para a superdifusão: 


y(t) = C Y1 

n=0 


/" n (Al)''/\"2' i 

r(no; + 1) 


OO 


71=0 


t na+ 1 (_l) n K 


2 n 


r(no; + 2) 


(1 < a < 2 


( 1 . 88 ) 


(1.89) 


'Uma função qualquer pode de ser escrita como uma série de Fourier, somatório de senos e 
cossenos e mais uma constante, desde que esta seja periódica, seccionalmente contínua e derivável 
exceto nos pontos de descontinuidade. Funções desse tipo podem ser utilizadas como condição 
inicial deste tipo de problema. 
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1.7 Fluxo Nêutrons X Cálculo Fracionário 

As equações cinéticas fracionais da difusão são uma abordagem útil para a descri¬ 
ção da dinâmica do transporte em sistemas complexos, que são regidos por processo 
anômalos da Difusão e por padrões de relaxamento não-exponencial. A difusão anô¬ 
mala pode ser modelada pela equação diferencial fracionária, tanto no tempo como 
no espaço. Para o espaço, parte do uso de divergência fracionária que modifica a 
expressão da difusão anômala, em a lei de Fick modificada . A aplicação desta di¬ 
vergência fracionária é apresentada em a definição do fluxo de nêutrons do reator 
nuclear. Quando se observa uma difusão anômala em escala de tempo, a modificação 
sugere o uso de equações cinéticas fracionadas. Os operadores de curvas fracionadas 
desempenharão talvez papel na teoria eletromagnética e Equações de Maxwell. [7] 

O perfil ou a cinética podem descrever melhor a realidade complexa. 

1.7.1 A Difusão de Nêutrons 

Os nêutrons movem-se em um reator em caminhos complicados em "zigzag" devido 
a colisões repetidas com o núcleo. Como conseqüência desse movimento, os nêutrons 
que estavam originalmente em uma parte de um reator e movendo-se em uma de¬ 
terminada direção com uma energia específica, aparecem mais tarde em outra parte 
do sistema, movendo-se em outra direção com alguma outra energia. Os nêutrons 
neste caso são ditos terem sido transportados da primeira região de energia para a 
segunda, e o estudo desse fenômeno é conhecido como teoria do transporte. 

A teoria do transporte é relativamente simples em princípio, e uma equação 
exata que governa os fenômenos de transporte pode ser facilmente derivada. Isso 
é chamado de Equação de Boltzmann, ou a equação de transporte, e o estudo da 
teoria do transporte é essencialmente o estudo desta equação. Infelizmente, é muito 
mais fácil derivar a equação de Boltzmann do que resolvê-la. Sob certas condições, 
no entanto, a equação simplifica consideravelmente e pode ser tratada de forma 
bastante direta. Esta versão simplificada da teoria do transporte é chamada de 
teoria da difusão e é o assunto do presente trabalho . Embora os requisitos para a 
validade da teoria da difusão raramente sejam plenamente realizados em problemas 
práticos do reator, o uso deste método normalmente fornece uma boa aproximação à 
solução de transporte exata e, por sua simplicidade, é comumente usado em muitos 
problemas de design de reator. [9] 

1.7.2 Conceito de Divergência Fracionária para Fluxo de 
Partículas 

Devido à relativa simplicidade e uso generalizado, a base da teoria local é dis¬ 
cutida primeiro. A teoria local faz uso de EDA (equação de difusão de advecção) 
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como: 


(1.90) 


dC 

~di 


\7.(—vC + DVC) 


onde C é a concentração de soluto, e D e v são tensores locais de dispersão e veloci¬ 
dade, respectivamente. O EDA baseia-se na definição clássica de divergência de um 
campo vetorial. A divergência é definida como a proporção do fluxo total através 
de uma superfície fechada para o volume fechado pela superfície quando o volume 
encolhe em direção à zero) como 


V.J = lim — [ J.ndS (1-91) 

v^o V Js 

onde J é vetor de fluxo, V é um volume arbitrário incluído na superfície Sené um 
vetor unitário normal a superfície. Isto é válido apenas se o fluxo for de fato uma 
quantidade representada por um vetor "ponto"em relação à escala de observação 
(por exemplo, fluxo de calor em material homogêneo). 

Então o limite existe e o operador reduz o produto de ponto familiar com o ve¬ 
tor gradiente [J-, -§:]■ A dispersão de soluções representa um contra-exemplo, 

devido principalmente às flutuações de velocidade que surgem apenas quando um 
espaço de observação aumenta, invalidando assim o limite. O fluxo de soluto ocorre 
devido aos efeitos combinados da velocidade média (advecção) e da flutuação da 
velocidade (dispersão). Os fluxos dispersivos para um determinado volume são tipi¬ 
camente calculados de alguma forma (volumétrico, estatístico) e aproximados pela 
primeira lei de Fick. Uma vez que a própria velocidade é uma função variável no 
espaço, à medida que o volume de controle diminui (como exige a divergência), as 
flutuações da velocidade e o fluxo dispersivo desaparecem. Portanto, a verdadeira 
divergência do fluxo de soluto macroscópico não poderá conter um termo dispersivo 
macroscópico. [9] 


1.7.3 Equaçao Cinética Fracionária 

As leis fundamentais da física são escritas como equações para a evolução do 
tempo de uma quantidade X(t) com dX ^ = — AX, onde esta poderia ser a equa¬ 
ção de Maxwell, a equação de Schrodinger, a lei de movimento de Newton ou as 
equações geodésicas. A solução matemática para o operador linear A é X(t ) = 
X(0) exp (—At), colocando r = A -1 a lei de relaxamento exponencial padrão é 

X{t) = X{0) exp (1.92) 

Sistemas complexos e investigações de suas propriedades estruturais e dinâ¬ 
micas foram estabelecidos na agenda da física. Essas estruturas com variações são 
caracterizadas através de: 
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ajiima grande diversidade de unidades elementares; 

b) fortes interações entre as unidades; e 

c) uma evolução não previsível ou anômala ao longo do tempo. 

Sistemas complexos e seu estudo desempenham um papel dominante nas ciências 
exatas e da vida, abrangendo uma riqueza de sistemas como óculos, cristais líquidos, 
polímeros, proteínas, biopolímeros ou mesmo ecossistemas. Em geral, a evolução 
temporal e dentro desses sistemas desvia-se das leis padrão correspondentes. Com 
o desenvolvimento de experimentos de maior resolução, esses desvios tornaram-se 
mais proeminentes. Pode-se ter um comportamento exponencial modelado como: 

X(í)=X(0)exp(-(^ ) (1.93) 

com 0 <a< 1, ou pode-se visualizar a lei do Poder assintótica como 

X(t) = vY(0) ^1 + , com n > 0 (1.94) 

Do mesmo modo, o processo de difusão em vários sistemas complexos ge¬ 
ralmente não segue as estatísticas gaussianas e, portanto, a segunda lei de Fick 
não descreve o comportamento de transporte relacionado. Especialmente observa-se 
desvios da dependência linear do deslocamento quadrático médio: 

<X 2 (t)> « lút (1.95) 

que é característico do movimento browniano, e uma consequência tão direta do 
teorema central do limite e da natureza markoviana do processo estocástico subja¬ 
cente. Em vez disso, a difusão anômala é encontrada em uma ampla diversidade de 
sistemas, sendo a sua marca registrada o crescimento não-linear do deslocamento 
quadrático médio ao longo do tempo, seguindo o Lei da Potência: 

<X 2 (Í)> « K a t a (1.96) 

o que é onipresente para diversos sistemas. Existe uma variedade de outros padrões, 
como a dependência do tempo logarítmico. O comportamento de difusão anômala 
manifestado em (1.91) está intimamente relacionado com a quebra do teorema do 
limite central, causado por distribuições amplas de cauda longa ou correlações de 
longo alcance. Essas amplas distribuições de saltos espaciais ou tempos de espera 
levam a uma forma de difusão não evolutiva do tempo, não gaussiana e possivelmente 
não-Markoviana, que se manifesta em fenômenos temporais não locais. Observe que 
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o coeficiente de unidade de difusão em (1.75) está tendo unidade [K a ] = cm 2 s~ a , 
de acordo com o expoente de difusão anômala. Esse expoente se 0 < a < 1 define 
o transporte sub-difusivo, e define fenômenos super-difusivos para 1 < a < 2. Para 
a = 1, o fenômeno de transporte é ordem inteira normal e Fickiano. 

A equação cinética de ordem inteira padrão quando integrada dá o seguinte 

Xi{t) - X 0 = -c f X(t) = -CoD^Xit) (1.97) 

J o 


D f _1 é o operador integral padrão de Riemann-Liouville. O número de densidade 
da espécie i, X t = X t (t) é uma função do tempo, e X t (t = 0) = X 0 é o número 
densidade de espécies i no tempo t = 0. Se retirarmos o índice i em (1.76) e substituir 
c por c a , então a solução da equação de difusão de ordem fracionária generalizada 
será: 

X(t) - X 0 = -c a í X(t) = —c~ a D^ 1 X{t) (1.98) 

J o 


m = 

k=0 


(—1 Yk)(ct) k ' a 
T(k.a + 1) 


(1.99) 


Pode ser escrita também de forma compacta, pelo uso da função Mittag-Leffler 
como: 


X(t) = Xo E a (-c a t a ) 


( 1 . 100 ) 


1.7.4 Descrição do Fluxo de Nêutrons em Reatores Nucleares 

O controle de um reator nuclear está intimamente ligado ao comportamento da 
população de nêutrons em seu núcleo, uma vez que são eles os responsáveis pelas 
fissões nucleares e, consequentemente, pela produção de energia. Sendo assim, é 
de fundamental importância o emprego de um modelo teórico capaz de descrever o 
comportamento físico desse sistema, como a Teoria de Difusão de Nêutrons. [8] 

Nessa teoria é feito o balanço da população em função do tempo, do espaço 
e da energia, considerando as possíveis fontes de geração e de perda de nêutrons. 

A descrição do balanço de neutrns no reator nuclear é definida pela teoria 
do transporte. As equações de transporte básicas são então aproximadas por várias 
equações diferenciais acopladas. Uma aproximação simplificada da representação 
do reator dada aos engenheiros é os conjuntos de equações de difusão de nêutrons 
em grupos de energia múltipla ou grupo de energia única. Em todas essas equações 
de difusão, o termo de vazamento tem a lei de difusão de Fick, onde o fluxo de 
nêutrons é assumido como uma quantidade pontual. Para uma maior representação 
do reator, várias dessas equações de difusão são formadas e modeladas por coefi¬ 
cientes de acoplamento região a região. A ciência da engenharia então prossegue 
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nestas aproximações para obter o modelo de função de transferência de reator e, em 
seguida, várias análises do sistema de controle são feitas. Para sistemas complexos, 
os modelos inteiros do reator podem não ser suficientes e, portanto, um modelo de 
ordem fracionada para obter perfil de fluxo ou cinética pode descrever melhor a 
realidade complexa. [7] 


1.7.5 Equaçao Clássica Constitucional da Difusão de Nêu¬ 
trons 


No sentido clássico, a equação constitutiva assume fluxo de nêutrons pontual, 
com v como a velocidade média dos nêutrons que atravessam uma área com n 
nêutrons por unidade volume como densidade de nêutrons. A quantidade de vetores 
que representam o fluxo de nêutrons é definida pela letra J. 

As equações seguintes irão elucidar as declarações clássicas. 


( 1 . 101 ) 
( 1 . 102 ) 

J= -D.V0 (1.103) 

Considere um volume fechado, a perda de nêutrons da superfície fechada é dada 
como a integral de superfície da corrente de nêutrons, J. dS. A perda que ocorre no 
volume por absorção é dada pela seção de absorção e, em seguida, tomando integral 
de volume de J2 a 4>dV. Essa perda total, quando equiparada ao termo fonte, dá a 
equação constitutiva clássica de difusão de neutrões, conforme descrito abaixo 


J = nv 
(j) = nv 


j s j dS+ j v j:, d v = j v 


= I SdV 
v 


A forma integral acima quando convertida em integral de volume é 

í (V.J + V <j>-S)dV = 0 
Jv a 

ou 

V. J + y (j) — S — 0, para o equilbrio 
Usando a expressão de J = —D.Vcf), obtemos o seguinte 

-D.V 2 0 + V 0-S = O 

‘ ^ n. 


(1.104) 


(1.105) 

(1.106) 


(1.107) 


estado estacionário, o (Right Hand Side ) RHS -lado direito- da equação cons¬ 
titutiva acima é zero, e se houver fluxo de mudança de tempo, então isso é colocado 
no RHS como 


D.W - V 


an 1 d(j) 
dt v dt 


(1.108) 
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onde: 

J2 a - é a seção transversal de absorção macroscópica [ cm _1 ]. 
n: nêutrons; 

S: é a força da fonte; 

cj): é o fluxo de nêutrons [neutron.cm~ 2 .s -1 ]; 

taxa de variação da densidade de nêutrons [nêutron.cm~ 3 .s~ 1 ] 

D: ê o coeficiente de difusão [cm] 

1.7.6 Discussão sobre Equações Constitucionais Clássicas 

A eqnação constitutiva de difusão de nêutrons clássica, conforme descrito, é 
baseada na divergência clássica da divergência de um campo vetorial. A divergência 
é definida como a proporção do fluxo total através de uma superfície fechada para 
o volume fechado pela superfície, quando o volume encolhe em direção a zero. 

V.J = lim — [ J.ndS (1.109) 

v^o 1/ J s 

Onde J é o vetor de fluxo, V é um volume arbitrário incluído na superfície S. O 
produto ponto do vetor J com a superfície dS é óbvio; isso é válido apenas se o 
fluxo é de fato uma quantidade de vetor "ponto"em relação à escala de observação. 
Difusão de nêutrons serve como contra-exemplo; Isto se deve principalmente, as 
flutuações de velocidade (mesmo em constante energia / temperatura) que surgem 
apenas quando o espaço de observação cresce, invalidando o limite. Também os 
nêutrons não estão mais em meio homogêneo. Os fluxos dispersivos para um 
determinado volume são normalmente calculados de alguma forma (volumétrica, 
estatística) e são aproximados pela primeira lei de Fick, como obtivemos ao derivar 

a equação constitutiva clássica para a difusão de neutrões, J = —D.Vcf). 

Um único nêutron em um meio de dispersão viaja em movimento em tons 
de ziguezague. Quando um grande número de nêutrons são considerados, existe um 
"movimento geral "da região de maior concentração (densidade) para a região de 
menor densidade (fluxo de números). Em casos hipotéticos, a velocidade permanece 
inalterada após a dispersão. Então, idealmente, em um meio de dispersão homogê¬ 
neo, o volume pode diminuir para zero (lim V —y 0), pois não haverá flutuações de 
velocidade. Apesar de considerar uma única fonte de energia, digamos nêutrons tér¬ 
micos, haverá flutuações de velocidade no ponto de medição, portanto, o fluxo não 
pode certamente ser denominado "quantidade pontual". À medida que o volume 
de controle diminui para zero, as flutuações de velocidade e o fluxo dispersivo de¬ 
saparecem. Portanto, em um verdadeiro sentido, o teorema da divergência clássica 
descere os efeitos reais da homogeneidade macroscópica e as flutuações associadas à 
difusão de neutrões em um reator. 
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Devido ao limite na definição de divergência, o teorema clássico de divergência 
de Gauss descontra o efeito de um grande volume até que o fluxo dispersivo possa 
ser aproximado por uma quantidade pontual. [9] 

A difusão ordinária começa a partir de duas equações. A primeira equação é 
a equação de continuidade: 


( 1 . 110 ) 


<90 d J 

d t dx 

Esta é uma lei de conservação fundamental relacionada com a variação tem¬ 
poral da densidade numérica ao fluxo macroscópico de partículas J. A equação cons¬ 
titutiva é a lei de Fick, uma equação fenomenológica que relaciona o fluxo com o 
gradiente de densidade (concentração) do número, através do coeficiente de Fick: 


J 



( 1 . 111 ) 


1.7.7 Divergência Fracionária em Equações de Difusão de 
Nêutrons 

dispersão de nêutrons através da superfície fechada, fechado por volume de 
controle finito, pode ser alterado para: 


—D 



d a 

dx a 



( 1 . 112 ) 


A expressão acima esta relacionada a primeira lei de Fick. Na qual [3 está entre 
1 e 2 e a está entre 0 e 1, números reais fracionários não inteiros. Nas formulações 
de divergência, o vetor 

d a d a d a 
dx a 1 dy a ’ dz a 

em que V" é o operador de divergência fracionária, que atua no vetor de corrente 
de nêutrons J. 



A equação constitutiva na forma de divergência fracionária, pode assim ser escrita: 


-D.V“</> + V 0-5 = 0 

^^ a 

(1.114) 

V“(-D.V0) + V 0-5 = 0 

(1.115) 

-D.V Q+1 0 + V 0-5 = 0 

^ J a 

(1.116) 


Convertendo para a forma unidimensional, para obter o perfil de fluxo de 
nêutrons em um reator, ou perfil de fluxo de nêutrons, perto de uma fonte em meio 
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de dispersão, analisaremos as soluções do formulário da seguinte forma: 


m+a: A 

b s^-£/+ s = 0 >-°<“< 1 (° 17 ) 

£•44 “V ^ + S = 0,...l</?<2 (1.118) 

d/X P ^CL 

Pode-se interpretar a forma simplificada de V o . J como um operador de diver¬ 
gência fracionária aplicado ao termo de dispersão de Ficki. 


1.7.8 Solução Clássica da Equaçao Constitucional de Difusão 
de Nêutrons para Ordem Inteira 

Este item servirá como uma revisão para a solução clássica simples da equação 
de difusão. Então, na próximo item, resolveremos a equação diferencial fracionada 
obtida. Isso demonstrará as variáveis espaciais em uma dimensão para simplificar. 


D. W - V 


d n 1 d(j) 

dt v dt 


(1.119) 


5 = 0 ( 1 . 120 ) 

D.V 2 <P - OCo - 1)^2 </> + s = ~ ( 1 . 121 ) 

a V CLZ 

O termo de fluxo é variável no espaço e no tempo; a constante k 0 0 representa o 
fator de multiplicação infinito[10]. A lei de multiplicação do termo original governará 
S. O método da separação de variáveis, dará a seguinte solução, para o termo de 
fluxo 

<t>(r,t) = <j>(r).<j>(t) (1.122) 


que pode ser substituído na equaçao constitutiva básica, e então surgirão as seguintes 
expressões: 

<j> = </>( r ).e ( " Aí) (1.123) 

D.V 2 (j){r) - (fcoo - 1)0(r) + S = — <p(r) (1.124) 

ZJ CL V 

A é positivo para estado sub-crítico, negativo para super crítico e zero para o reator 
de equilíbrio crítico. Substituímos a coordenada espacial r, por x e com a subs¬ 
tituição de B por encurvamento geométrico, seguimos a forma simples. A solução 
temporal é evitada por simplicidade. 


B 2 


(goo ~ JQZla + f 
D 


(1.125) 
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d 2 (f>(x ) 
dx 2 


+ B 2 .(f){x ) = 0 


( 1 . 126 ) 


O termo 5 2 é conhecido por flambagem geométrica, pois este termo se relaciona 
com curvatura geométrica do perfil de fluxo de nêutrons no espaço: 


B 2 = 


d 2 cf)(x) 
dx 2 

4>{x) 


( 1 . 127 ) 


quando este é igual à deformação de material, o reator está em equilíbrio crítico. 

Aqui podemos aplicar o método de Laplace padrão com condições de con¬ 
torno iniciais em x — 0 no ponto central da geometria do reator com fluxo constante, 
e nas paredes em x = a com fluxo zero. A fórmula geral após aplicação da transfor¬ 
mada de Laplace para a derivada da função é indicada abaixo e é aplicada para se 
ter uma forma polinomial 

C{<p(x)} = $(s) 

i 


k =0 


s 2 $(s) - 

s 2< f>(s) — 


.d 2 k 1 <f>(x) 


dx 2 ~ k ~ 1 


+ 5 2 $(s) = 0 


( 1 . 128 ) 

( 1 . 129 ) 


x=0 


Como as condições iniciais 


d(f>(x) 

— S(f>(x ) 

x=0 

dx 

d(f)(x) 

dx 

= 0 e 4>(x) 


+ B 2 $(s) =0 


( 1 . 130 ) 


a ?=0 




= C, temos: 


x=0 


De onde obtemos 


s 2 $(s) + sC + B 2 $(s) = 0 

, sC 
$(s) = 


+ B 2 


Aplicando a transformada inversa de Laplace, tem se: 

4>(x) = C cos (Bx) 


( 1 . 131 ) 

( 1 . 132 ) 

( 1 . 133 ) 


A solução para a equação de ordem de número inteiro clássica também é obtida 
da seguinte maneira: 


d 2 (p( 


x) 


dx 2 


+ B 2 .cf)(x ) = 0 


Considerando a simetria: 


4>(x) = d\ cos (Bx) + a 2 sin (Bx) 

4>{x) = (f>(—x) 

4>(x) = ai cos (Bx) 


( 1 . 134 ) 

( 1 . 135 ) 

( 1 . 136 ) 

( 1 . 137 ) 
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Para avaliar B, a condição de fronteira (limite) para x = a (onde a distância é 
extrapolada) implica qne o fluxo torna-se nulo. Assim 

0(±o) = ü\ cos (Ba) = 0. (1.138) 


Como Oi nao pode ser zero, cos (Ba) = 0, segue que B pode ter valores iguais a 


r> _ mr 

r>n ~ 2a ‘ 

Assim 


fmrx\ 


(x) = C. cos J , (n = 1,3,...) 


(1.139) 


B 


2 

1 


' 7T 1 2 

.2 a. 


(1.140) 


Os outros termos de ordem superior são chamados de autovalores e seu termo de 
fluxo relacionado decai com o tempo; apenas o primeiro termo dominante governa 
o perfil de fluxo. 


1.7.9 Solução da Divergência Fracionária baseada na Equa¬ 
ção de Difusão de Nêutrons ( para Ordem Fracionária) 

Extendo o método acima, resolveremos a equação diferencial fracionária: 


d /3 (f>(x) 

dxP 


+ B 2 4>(x ) 


0, 


1 < p < 2 


(1.141) 


s^(s 


d a f(x ) 
dx a 

d^~ l cf)(x 


n— 1 


} = s«F(s) - ^ 


,d 


a—k—1 


f(x) 


k =0 


x=0 


dxP~ l 


dP~ 2 (j)( 


X) 


dxP~ 2 


x=0 


d a ~ k - 1 
+ 5 2 <í>(s) = 0 


x=0 


(1.142) 

(1.143) 


A equação (1.142) expressa acima, trata-se da transformação de Laplaee para a 
definição LHD da derivada de ordem fracionária. Nesta expressão, o segundo e o 
terceiro termo do lado esquerdo têm derivada fracionária do fluxo no ponto inicial, 
que é fisicamente difícil de definir, e realizar o mesmo por medidas experimentais é 
difícil nesta fase. Vamos tentar fazer uso da transformação de Laplaee da definição 
do RHD de Caputo, conforme indicado abaixo: 




d a f(x ) 
dx a 


n— 1 


}=ü“F(s)-^ 


k =0 


s^(s) - sP-^iX) 


— S ' 


0_ 2 d<j)(x) 


X=0 


dx 


x=0 

(1.144) 

+ 5 2 $(s) = 0 

(1.145) 


x=0 
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Nós relacionamos a expressão acima fisicamente com a condição inicial anterior 
tomando o segundo termo como C e terceiro termo como zero como feito no item 
anterior. Aqui, a derivada de ordem inteira vem como condição inicial, portanto, 
fisicamente realizável a partir de medidas e observações. 

s P~ l C 

*(*) = (!- 14 6) 


(j>{x) = ccr 1 


' s^C ' 
sP + B 2 


(1.147) 


A solução da equação diferencial fracionária para a equaçao constitutiva de equilíbrio 
de nêutrons, é expressa da seguinte forma, usando a transformada de Lapalace: 


De onde obtemos 


C{{E a (-\t a )} 


s° + A 


^(x)=^(-5V)=C^ 

k =0 


(-sV) fc 

fpDl) 


(1.148) 


(1.149) 


<j>{x) =C + C 


-bV) 


+ C 


■bV) ! 


+ c 


-B 2 x py 


+ c 


-B 2 x 0 y 


r(/3 +1) r(2^ + i) r(3 / d + i) r(40 + i) 


+ 


(1.150) 
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As equações representadas acima representam o mapeamento de fluxo ob¬ 
tido pela solução de equações constitutivas de nêutrons de ordem fracionária, que 
são obtidas pelo conceito de divergência fracionária. O gráfico da figura abaixo re¬ 
presenta, o comportamento para fluxo quando variamos o valor de x e de k para 
13 = 1,5 ;/3 = 1,8 e $ = 2,0. 



x 


Figura 1.3: Mapeamento do fluxo de nêutrons onde C(x) = 4>(x) 


Observemos o padrão de fluxo clássico e o padrão de fluxo de ordem fracionária 
quando tomamos a ordem fracionária j3 igual a 2. A solução na forma clássica 
relacionada com a função cosseno tem sua representação em série como: 


o 4 fi 


cos(*) = l-- + 


P = 2 , 


3=2 (x) = CE i3 (-B 2 x 2 ) = CJ2 


(~B 2 x 


2 ^ 2 \k 


=2ÍX) = C 


r <2fc + 1) 

j (B 2 x 2 ) ( B 2 x 2 f (. B 2 x 2 f 


r(3) r(5) 


r(7) 


como T(n + 1) = rd, temos: 


(1.151) 


(1.152) 


(1.153) 


4>/3=2(x) = C 


{.B 2 x 2 ) 
2 ! 


+ 


(BV 2 ) 

4! 


(. B 2 x 2 ) 3 


~ C cos (Bx) 


(1.154) 
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Portanto, quando a ordem fracionária é igual a ordem inteira, obtemos o perfil de 
fluxo clássico. Isto é prova de nossa suposição de que, de fato, o fluxo de nêutrons, 
não sendo uma quantidade pontual, pode ser representado como uma divergência 
fracionária de uma ordem menor que a unidade. Srivastava tentou ter funções 
trigonométricas realizadas pela função Mittag Leffler e usou as variantes da função 
Mittag-Leffler para resolver equações diferenciais fracionárias. O cosseno fracionário 
e as funções do seno fracionário são manifestações dessa importante função (Mittag- 
Leffler) no cálculo fracionário. Na solução analítica acima, demonstramos o uso 
de RHD para derivativos fracionários para nossa conveniência; No entanto, o LHD 
também é usado quando as condições iniciais das derivadas fracionárias são bem 
compreendidas e podem estar relacionadas à física. Atualmente, nossas ferramentas 
matemáticas vão muito além dos limites de nossa compreensão física. Assin sendo, 
qual a definição - RHD ou LHD - aplicável depende dos limites de nossa compreensão 
física. Por exemplo, nos estudos de polímeros, os modelos viscoelásticos tomam LHD 
e relacionam os derivados fracionários da variável de tensão com o estresse e obtêm 
soluções para o relaxamento do estresse. A extensão da pesquisa é necessária para 
obter informações sobre as condições iniciais baseadas em derivadas fracionadas e 
para atribuir significância física, quando a LHD é usada para solução; Atualmente, 
a teoria matemática vai muito além do nosso entendimento físico. 
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A Figura (1.3) dá vários valores para a ordem fracionada f3. Os valores da 
função Gamma são obtidos a partir de uma calculadora padrão da função Gamma. 
Os valores utilizados são T(2.5) = 1.32934, T(5.5) = 52.3427 r(2.8) = 1.6765, 
r(4.6) = 13.3813, r(6.4) = 240.8337, T(8.2) = 7562 


X 

0 = 1.5 

CO 

II 

02 . 

0 = 2.0 

0,0000 

1,0000 

1,0000 

1,0000 

0,3996 

0,8257 

0,9214 

0,9238 

0,7853 

0,5510 

0,6442 

0,7071 

1,1780 

0,2744 

0,3242 

0,3820 

1,5700 

0,0000 

0,0000 

0,0000 


Figura 1.4: Solução para equaçao constitutiva fracionária para várias ordens. 

Por simplicidade e apenas para representação, B = 1, C = 1 são tomados 
e apenas quatro termos na série de soluções são considerados. O ponto x = 0 é o 
centro do reator e x = 1,57 é o limite extrapolado do reator nuclear no qual o fluxo 
é zero. A tabela mostra o fato de que quando a ordem fracionária tende para 2, a 
equação constitutiva dá uma distribuição do cosseno do perfil de fluxo. 


1.7.10 Difusão de Nêutrons em Medição Nao-Multiplicativa 
em um reator. 


Pense numa fonte de nêutrons no plano mantida ao longo do eixo Y. A situação 
é semelhante a considerar um anteparo mantido na borda do reator (um meio de 
multiplicação). Anteriormente , o exemplo foi a difusão de nêutrons na multiplicação 
de média de -a"para "+ a", e no reator aberto a distribuição de fluxo foi obtida. 
Asumiremos o fluxo de nêutrons finito em x = 0 e zero em outro lugar. Também 
em x = 0, o fluxo de nêutrons tem um gradiente finito (não nulo). A equação 
constitutiva de ordem inteira clássica, é expressa por: 


D 


torna-se: 


d 2 (p(x) 
dx 2 

d 2 cj)(x 


E *(*)+ s 

< n. 


1 d(j) 
v dt 


D 


dx 2 


^ = 0 

< * a 


(1.155) 


(1.156) 


na falta de um termo de multiplicação de fonte; o último pode ser simplificado como: 


D d 2 (j)(x) 
dx 2 


k 2 cj)(x ) = 0 


(1.157) 
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Usando a Transformada de Laplace, temos: 


s 2 <f>(s) 


d(j)(x) 


dx 


x=0 


S<t>(x ) 


+ k 2 $(s) = 0 

x=0 


Considerando: 


d<f)(x) 

dx 

1? 

co 

• 1 

= — Ck ; 4>(x) 


O 

II 

O 

II 


Obtemos: 


sC^kC _ sC kC 
^ ~~ s 2 - k 2 ~ s 2 - k 2 ~ s 2 - k 2 
Aplicando a transformada inversa de Laplace, obtemos: 


4>(x) = C cosh (kx) + C sinh (kx) = Ce kx 


(1.158) 


(1.159) 


(1.160) 


A solução é mostrada propositadamente para ter termos hiperbólicos. A solu¬ 
ção com frações parciais teria dado a decomposição exponencial direta em relação 
à x. Isso é propositalmente mostrado para explicar o seguinte com as expressões 
constitutivas de divergência fracionada. A equação não local e suas transforma¬ 
ções de Laplace na definição de Caputo, conforme feito anteriormente, produzem os 
seguintes resultados: 


p d p <i>{x) 

dx^ 


k 2 (p(x) = 0 


(1.161) 


s^(s) 


d(j)(x) 


dx 


x=0 


S(j)(x) 


+ k 2 $(s) 

x=0 


Considerando: 


d(f)(í r) 
dx 

1? 

co 

j 1 

= —Ck] 4>(x) 


O 

II 

O 

II 


Obtemos: 


0 




sC — kC 

S 2 — k 2 


sC _ kC 

S P _ ]ç2 gp _ fc 2 


(1.162) 


(1.163) 


4>(x) = CR/ 3 ^[k 2 , x] + kCFp[k 2 , x 


(1.164) 
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Expandindo a forma compacta das funções ReF em notações em série, obtemos: 


<Hx) = cJ2 

k =0 


rjç2\n x (n+l)l3—l 

r({n + l}/3-l) 


(1.165) 


A equação clássica com ordem inteira foi resolvida mantendo a soma das funções 
hiperbólicas. Ambas as funções divergem e são ilimitadas para grandes distâncias da 
origem, mas o resultado final produziu um resultado convergente. Isto é exatamente 
o que é enfatisado para as duas soluções de série divergentes na forma de funções F 
e R, para a equação diferencial de ordem fracionada. A soma líquida deve chegar 
ao resultado da convergência. Essas funções generalizadas são traduzidas para ex¬ 
pressões hiperbólicas, trigonométricas e outras expressões exponenciais para índices 
adequados. 


Conclusão 

Aqui, assinalamos o uso da definição de divergência fracionária. Também foi 
discutido um método não local de obtenção da quantidade de fluxo de nêutrons 
a ser aplicado, pois argumenta-se que o fluxo de nêutrons de ponto-quantidade 
clássico não reflete a realidade. O fluxo de nêutrons não é uma quantidade pontual 
e, portanto, a definição de divergência fracionária dá a imagem real. Além disso, 
não só o método não local de obtenção de divergência dá efeito às flutuações da 
velocidade, mas também capta o efeito de proximidade real da presença de nêutrons 
vizinhos. O conceito aqui apresentado requer experimentação e técnicas de regressão 
numérica para mapear o fluxo de nêutrons na ordem fracionada. Além disso, a ordem 
fracionada será diferente de reator para reator devido à heterogeneidade, variando 
com as configurações do reator. Este método de perfil de fluxo visará o controle 
eficiente das usinas nucleares. 
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